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ABSTRACT 


This  report  investigates  the  difficulties 
encountered  with  respect  to  a certain  linear  least 
squares  fit  problem  which  arises  in  the  field  of 
holographic  interferometry.  It  was  found  that  the 
addition  of  another  variable  to  the  equations  of 
this  linear  least  squares  fit  problem  explained 
these  difficulties  and  provided  consistent  results. 
This  report  also  demonstrates  the  usefulness  of 
Stewart's  recently  published  perturbation  bounds  for 
the  purpose  of  error  bounding. 


ADMINISTRATIVE  INFORMATION 


This  work  was  performed  under  NAVSEA  Mathematical  Sciences  Program 
"Numerical  Methods  for  Naval  Vehicles,"  Program  Element  61153N,  Task 
Area  SR  0140301,  Task  13321,  DTNSRDC  Work  Unit  1-1808-010.  The  NAVSEA 


cognizant  program  manager  is  Ms.  B.  Orleans,  SEA  03R 


INTRODUCTION 


Recently  Dr.  Surendra  K.  Dhi r (Head,  Numerical  Structural  Mechanics 

Branch,  Code  1844)  called  the  author's  attention  to  difficulties  arising 

from  a linear  lea»L  squares  fit  problem  involving  holographic  displace- 
1* 

ment  measurements.  In  the  problem  described,  it  is  desired  to  fit  a 
polynomial  of  the  form  w = ax  + by  + cz  by  a set  of  experimentally  ob- 
tained data  points  (w. , a^,  b^,  c^),  i=l,...,n;  but  under  certain 
circumstances  the  polynomial  would  have  to  be  set  up  in  a differenced 
form,  i .e . , 


Although  the  polynomial  can  be  fitted  to  this  "differenced"  data,  the 
relationship  between  this  "differenced”  problem  and  the  original  problem 
is  neither  apparent  nor  simple. 

One  possible  approach  is  to  regard  the  two  least  squares  fit  problems 
as  perturbations  of  one  another  (after  adding  a redundant  equation  to  the 


"differenced"  system).  It  was  with  this  approach  in  mind  that  Dr. 

Elizabeth  H.  Cuthill  (Assistant  for  Numerical  Analysis,  Code  1805) 

suggested  investigation  of  Stewart's  recently  published  perturbation 

2 

bounds  for  the  linear  least  squares  fit  problem.  These  bounds  were 
found  to  be  very  useful  indeed.  This  report  documents  their  application 
to  the  problem  in  holographic  interferometry. 

STEWART'S  PERTURBATION  BOUNDS 

Let  A be  an  m x n real  matrix  of  rank  r.  It  is  possible  to  find 
orthogonal  matrices  U and  V of  orders  m and  n,  respectively,  such  that 


1 


1 


I 0 


T 

U AV  = 


where  i 


0.  (it  will  be  remembered  that  an  orthogonal  matrix 
T 


II  is  a real  matrix  whose  transpose  U is  its  inverse.)  The  product 
T 

U AV  is  referred  Lu  as  the  singular  value  decomposition  reduced  or 


standard 
A.  The  product 


form  of  A.  The  a.,  i=l 


,r  are  called  the  singular  values  of 


1/a 


1 /<J, 


l/<) 


is  called  the  pseudo- inverse  of  A,  written  A . 


he  pseudo-inverse  has 


many  of  the  properties;  of  the  inverse  of  a square  matrix.  It  is,  in 
fact,  a generalization  of  the  concept  of  a matrix  inverse.  For  more 
information  on  these  topics,  the  reader  is  referred  to  Stewart. 

The  formulas  for  Stewart's  perturbation  bounds  require  that  A has 
been  reduced  to  the  singular  value  decompos i t ion  standard  form.  This 


t 


reduction  involves  no  Loss  of  generality  since  any  m x n linear  least 
square  fit  problem  Ax  = b and  any  of  its  perturbations  (A+E) (x+h)  = b+k  can 
be  reduced  to  this  required  form  by  the  respective  transformations 


(UTAV)(VTx)  » 
(llTAV+UTEV)  (Vrx+V*h) 


T 

Ob 


T T 
U b+U  k 


For  the  purpose  ot  our  problem  we  shall  also  assume  that  the  columns  of  A 
are  linearly  independent  and  that  the  first  n rows  of  A are  a^ so  linearly 
independent . 

Let  the  linear  least  squares  lit  problem  A(x+h)=b+k  result  from  a 
perturbation  of  the  right-hand  side  of  Ax  = b.  Define  < as  ||  A 1 1 |!  A ||2 

and  n as  | bj ':|  , / 1|  A|;  , , x | .,  where  bf  is  the  sub-vector  of  the  first  n 
components  of  b.  Define  kj  similarly.  Then  Stewart  provides  the  pertur- 
bation bound 

h ||  2/  ||  x,  , - , kj  i,2  / b1  i2  (Inequality  1) 


Next  let  tu.  linear  least  squares  lit  problem  (A+E)x  = b result  from 
a perturbation  of  the  coefficient  matrix  of  Ax  * b . Define  k as 
||  A ||  0 ||  , where  B^j  is  the  inverse  of  the  nt'1  order  principal  sub- 

matrix ot  A+E.  Let  E be  the  n1*1  order  principal  submatrix  of  E and  let 
E,j  be  the  (m-n)  xn  submatrix  ol  E below  E j ^ . Define  n and  bj  as  before. 


The  reader  is  no  doubt  familiar  with  the  concept  of  the  norm  of  a 
vector  x,  written  ||x  j[,  and  defined  by 


r n 

r E 
i=l 


where  the  x.  are  the  components  of  the  vector  x.  Since  there  are  several 
norms  defined  for  vectors,  let  us  be  more  precise  and  refer  to  the  above 
norm  as  Ijx^.  For  A,  an  m x n real  matrix,  ||  A ||  2 i-s  defined  to  be  o^, 

the  largest  singular  value  of  A.  It  can  be  shown  that  1 1 x 1 1 2 = i | V x 1 1 2 and 
||  A ,=  II  UT  A V ||  2 where  U and  V are  unitary  matrices.  Accordingly  this 

norm  (along  with  others)  is  said  to  be  unitarily  invariant.  For  more 
information  on  norms  and  their  properties  the  reader  is  referred  to 
Stewart . 


I i 


Let  b^  be  the  m-n  subvector  of  b below  b^.  Then  Stewart  obtains  the 
perturbation  bound 

llhl'  _ 


(Inequality  2) 


Define  the  real  valued  matrix  function 

i^(F)  = — 


(i+  ]|f||; 


Stewart  bounds  the  last  term  of  Inequality  2 and  obtains  another  pertur- 
bation bound 


< K 


lEllll2  Mb2 
+ k n 


INI- 


?!l  J 

ill,  \ l|A||  / 


( Inequai i ty  3) 


Clearly  these  three  inequalities  can  be  used  to  provide  "worst  possible 
case"  bounds  for  a Riven  least  square  fit  problem  Ax=b.  In  the  case  of 
Inequality  1 one  needs  only  to  provide  some  upper  bound  for  ||k  ||2  to 
obtain  such  a bound  for  |j  h | j.,  / ||  x ||  when  only  the  right-hand  side  b is 
perturbed.  As  regards  Inequalities  2 and  3 one  must  provide  upper  bounds 
f°r  'lEii^l  and  ||E?1 1|  to  obtain  similar  worst  possible  case  bounds  for 
/!|*!l,  when  only  A is  perturbed.  The  reader  is  reminded  that 


II  (A 


n + Eii> 


S l|An"2 
l - II  a" 


11 


l2  He 


(see  Wilkinson,4  p.  92)  if  ||A~|  ||  , II  Eu  II  2 


112 

< 1. 


In  Inequality  3 the 


function  ip  may  be  bounded  from  above  by  1.0. 


/ 


THE  TESTS 


We  shall  refer  to  the  7x3  linear  least  squares  fit  problem  Ax  = b 

where 

- . 21458051E+00  - . 59223911E+00  - . 18537080E+01 

- . 35869262E+00  - . 54449323E+00  - . 18176699E+01 

- . 22761852E+00  - . 47379546E+00  - . 18959881E+01 

A = - . 30491734E+00  - . 45507956E+00  - . 18743727E+01 

- . 37405904E+00  - . 43714104E+00  - . 18487759E+01 

- . 23612950E+00  - . 35108814E+00  - . 19235878E+01 

- . 38382191E+00  - . 32889374E+00  - . 18685387E+01 

and 

. 36055440E-03 
. 35650322E-03 
. 35852881E-03 
b = . 35650322E-03 

. 35447764E-03 
. 35650322E-03 
. 35245205E-03 

as  the  "original"  problem.  It  is  estimated  that  the  elements  of  the 
columns  of  A and  b are  subject  to  uncertainties  of  5%,  11',  2%,  and  5%, 
respectively. 

We  define  the  "differenced"  7x3  linear  least  squares  fit  problem 
thus:  for  i = 1,...,6  subtract  the  i1*1  equation  from  the  i+lst  equation 
of  the  original  problem;  the  seventh  equation  of  the  differenced  problem 
is  the  difference  of  the  first  and  seventh  equations  of  the  original 
problem.  We  define  the  "augmented"  7x4  linear  least  squares  fit  problem 
merely  by  adding  a fourth  column  of  1.0's  to  the  coefficient  matrix  of 
the  original  problem  and  another  unknown  x^  to  the  solution  vector.  Now 
solve  this  augmented  problem  for  the  value  of  x^ , take  the  original  7x3 
linear  least  squares  fit  problem  Ax = b and  subtract  the  value  of  x^  from 
the  right-hand  side  b of  the  original  problem.  We  shall  refer  to  the 
resulting  7x3  linear  least  squares  fit  problem  as  the  "modified"  problem. 

These  linear  least  square  fit  problems  were  perturbed  ten  times  each. 
The  ratio  of  the  norm  of  the  solution  perturbation  to  the  solution  norm 
was  computed  along  with  their  bounds  provided  by  Inequalities  1 through  3. 
These  results  are  compared  in  Tables  1-5  for  perturbations*  of  either  or 
* The  perturbations  in  Table  5 were  reduced  by  a factor  of  10  to  meet  the 
condition  ||A^|||  ||e^||<1  for  the  approximation  of  ||(A^  + E ) ‘ jj  } . 
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mm 


TABLE  5 - UPPER  BOUNDS  TO  INEQUALITIES  2 AND  3 BOUNDS 


both  the  right-  and  left-hand  sides  of  these  problems.  We  shall  comment 
on  these  results  in  the  next  section.  These  results  were  obtained  by 
single  precision  computation  on  the  CDC  6400  at  DTNSRDC.  Subroutines 
from  the  I MSI.  Library"*  were  used  to  compute  singular  value  decompositions, 
solve  linear  least  squares  fit  problems,  and  generate  random  numbers. 


21458043E+00  - . 6161 7945E+00  - . 18521478E+01 
36091016E+00  - . 52242982E+00  - . 1 8192 5 39E+01 
21953501E+00  -. 52034 14 1 E+00  - . 18609241E+01 
31 149031E+00  -.415637 15 E+00  -. 1880861 1E+01 
36469258E+00  - . 4 3905926E+00  - . 18250381E+01 
23855953E+00  -. 324827 39 E+00  -. 19453948E+01 
3829  7 324E+00  - . 33146348E+00  - . 18640814E+01 


. 36055496E-03 
. 34  7 1 2449E-0  } 
. 3589  3074E-0 3 
. 34162634E-03 
. 35680069E-0 3 
. 34089 129E-03 
. 35576825E-0  1 


be  the  perturbations  oi  the  A and  b matrices  defined  previously.  The 
solution  of  the  original  problem  Ax = b is 


- . 41 300149E-04 
-.66551 182E-04 
-.  1679  30  151'.- 03 


The  solution  x'  ot  the  problem  Ax'  = b (only  b perturbed)  is 


-.677  666  79K-04 
-.771620851-04 
-. 1 585 1 240E-0 3 


The  solution  x"  of  A'x"  =b  (only  coefficient  matrix  perturbed)  is 


-.  56582677 1.-04 
- . 72542926E-04 
-. 1 6455871 E-0 3 


The  solution  x' 


b ' (both  b and  coefficient  matrix  perturbed)  is 


-.8812058 3E-04 
100616.  >K-0  1 
-.  1 50074061.-0  1 


These  results  should  prove  useful  us  benchmark  results  for  those  who  wish 
to  run  the  problems  published  in  this  report. 


!: 


h 


- 
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OBSERVATIONS  AND  CONCLUSIONS 

Let  us  begin  by  explaining  the  discrepancy  between  the  solutions  of 
the  original  and  the  differenced  linear  least  squares  fit  problems.  A 
system  of  m linear  equations  in  n unknowns  Ax  = b is  said  to  be  consistent 
if  there  exists  a solution  vector  x which  satisfies  each  equation  of  this 
system  exactly.  It  is  obvious  from  the  rules  of  algebraic  manipulation 
that  any  solution  of  a consistent  system  is  likewise  a solution  of  its 
differenced  system.  If  Ax = b is  inconsistent,  we  define  the  linear  least 
squares  fit  "solution"  to  be  that  real  vector  x which  minimizes  the  norm 
of  the  residual  vector  r = Ax  - b . Thus  it  is  merely  a notational  con- 
venience to  write  an  inconsistent  linear  least  squares  fit  problem  as 
Ax  = b,  since  equality  obtains  only  when  we  include  the  residual  vector  r 
thus 

Ax  = b + r 

If  we  construct  the  differenced  system  by  subtracting  the  equations  of 
this  last  system  from  one  another,  we  can  see  that  the  residual  vector  r' 
for  x with  respect  to  the  differenced  system  is 


In  general  r'  is  not  the  residual  vector  of  minimum  norm  for  the 
differenced  system  and  hence  x is  not  the  least  squares  fit  solution  for 
the  differenced  system. 

Next  let  us  consider  the  modeling  of  the  original  7x3  problem.  If 
we  include  an  extra  variable  k in  each  equation  of  the  original  system 
(thus  obtaining  the  7 x A augmented  system),  we  find  that  the  first  three 
unknowns  of  the  augmented  system  least  squares  fit  solution  do  not  agree 
at  all  with  the  least  squares  fit  solution  of  the  original  system.  Rather, 
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1 

they  coincide  with  the  least  squares  fit  solution  of  the  7x3  system 
whose  coefficient  matrix  is  A but  whose  right-hand  side  is  diminished  by 
the  value  of  k (modified  system).  Clearly  this  fourth  unknown  k is 
crucial . 

Finally,  this  investigation  has  demonstrated  the  value  of  using 
Stewart's  perturbation  bounds  in  the  error  analysis  of  the  linear  least 
squares  fit  problem.  In  particular  the  Inequality  1 estimates  seem  to 
bound  ||  h |!  ,/ |j  x || ,,  quite  nicely  when  only  b is  perturbed.  In  fact  these 
Inequality  1 bounds  are  roughly  a third  of  the  worst  possible  error 
estimate  afforded  by  Inequality  1.  As  might  be  expected  the  Inequalities 
2 and  3 estimates  do  not  bound  ||h||_/||x||  l when  A is  perturbed  as  closely 
as  the  Inequality  1 estimates  do  for  perturbed  b.  However  at  worst  the 
Inequalities  2 and  3 estimates  appear  to  be  not  greater  than  eight  times 
||h|L/||x||  . Usually  the  Inequality  3 bounds  are  "sharper"  (closer)  than 
the  Inequality  2 bounds,  although  this  is  not  always  the  case.  These 
Inequalities  2 and  3 bounds  are  roughly  a tenth  of  the  worst  possible 
error  estimates  provided  by  Inequalities  2 and  i,  respectively. 

So  far  all  results  obtained  indicate  that  the  errors  from  a joint 
perturbation  of  A and  b are  roughly  the  sum  oi  the  respective  errors  and 
hence  may  be  hounded  by  the  sum  of  the  bounds  provided  bv  Inequalities  1 
and  2 or  1 and  3. 

P ROC  RAMS 

The  following  programs  were  used  to  obtain  some  of  the  results 
published  in  this  report.  The  reader  is  reminded  of  the  special  row  and 
column  rank  conditions  holding  for  this  problem.  In  particular,  the 
procedure  for  obtaining  an  mxm  unitary  matrix  U'  from  the  mxm  unitary 
matrix  U output  from  the  I.SVALR  subroutine  works  only  under  these  condi- 
tions. (The  last  m-n  columns  of  the  mxm  identity  matrix  1 are  annexed 
to  the  mxm  matrix  U.  U'  is  obtained  by  orthogonal  izing  these  m columns 
using  Gratn-Schmidt  algorithm.  Of  course,  this  procedure  can  be  generalized 
by  modifying  it  with  suitable  before  and  after  permutations.) 


P7D67AM  SKDP703  OUTPUT,  TAPE6=0UTPUTI 


DIMENSION  A (7),  Dll ) ,C  « 7 1 , A A < 7,  3 I , BB(  71  , * ( 7 I , N<  ARE  A < 99)  , X ( 7 I , Y < 7»  * 
L 7171  , X N (7  I ,S(7),CC(7),Cu0S(3«3),Q(7,3l,H(7l«tf  (7,7)  » & ( 7 , 3) 
DIMENSION  TK<(3.  3 I, B0(3> 

DIMENSION  T (7  I 

DIMENSION  QQI7.3I  .38(71 ,0C( 7,7> ,aD (3,3 ».a£ <31 
DIMENSION  QUI7.7I  ,X KK<7 I , XK ( 71 
DIMENSION  OP  I 71 
DIMENSION  QPI7.3I 

DIMENSION  Ell  (7,3  I,  G&(7,3l  , GO  (7),  GO  (7, 31  , 3Y  (3,  31  ,G4  (7 ,3)  ,ELC  7,  2>  , 
l 07(7, 31,  ET  (7,31 
DIMENSION  E 21 (4, 3 I 
DIMENSION  CFI7I 

00  3%  JO-1, 10 

IP  ( JO.EQ. 1)  »U=JQ 

XX=-S. 165*6.1 
Yf=-2. 1*59-0  .7 
22=2.53-19. 4 
XO=VO=23=0 .0 

R9=Sa7T((XX-XOI*(XX-XOI*(YV-YOI*( YY-VGI *CZZ-Z0l*<ZZ-ZU)l 

X(l»=X( 3>  = X (61  = 4.25 

X(2I=X(5I=X(7|=5. 7 

X (.1 =5.-6 

Y(1I=Y(2I*1.87» 

Y(JI=Y(4I=Y  (5l  = ).25 

Y(6I*Y<71=-1.25 

Z ( 1 > =7  (3I  = Z (6M-5.12 

7(21=2  I5)=Z  (71*-5.2 

Z(4l=-5.63 

00  l <=1,7 

7(<»=Sa7T((XX-X(<n*(KX-X(K»)*(YY-Y(<n*(YV-Y(<M*(ZZ-Z(K»)  = 

L ( ZZ-2 ( < 1 1 ) 

A(<I=(XX-X0I/RJ*(XX-K<KII/R(KI 
9(KI  = <YY-Y(,I/RJMYY-Y(KM7R(KI 
1 CI<I=(7Z-Z0I/RJ»(2Z-2<KII/R(KI 
XN( II =17.8 

XN( 21 *XN (4 l =XN( 51  =1 7.6 
XX ( 31 =17.7 
XN( il =17.5 
X N ( 71 =17.4 

XlAM33A=5145.0*3.937  E-09 

00  *7  1=1,7 

oa(i.ii=A( n 

00(1,21=6(1) 

aai i, 3i=ci i i 

09(11 = X L AX  80A  * X X ( 1) 

Cr (II =38(1  I 

47  OP(H  =09  (I  1 

9XR1=SaPT(O8(ll*»2*QBI2)**2»OB(3l**2*09(4l*=2*QB(5»**2» 

1 33(51 ••2*03(71  **21 
00  *8  J* 1, 3 
DO  41  1=1, 7 
02(1, JI=Q0(I, Jl 

48  aj(  I,  Jl  =Q0  ( I,  Jl 


n = 7 
N*  3 
ISi«l 

MRITE(b,99l ((Q3(I,J>,  J=1.3I  ,1*1,71 
99  F3R1AT(lX,2MaQ,  3E16.8I 

CALL  l5VAlR(QQ,1, N, M,  N, I$W,WKAREA,QE,aC»3DI 
A9RH  = A MA  XI  (Q£(l),QE  (21  ,4£(3)  I 
XA=»iI2  = AMAXK1.17QE (U .l./OE (21  , 1./QEI3H 
MRITE(b,62l ((QC(I,JI, J=l,3> ,1=1,71 

62  F0R«AT(lX,2naC,3Elb.7l 
CALL  tfCOMPCQC.aJ.N.NI 

WRITE  16,63)  ( ( QJ  ( I »JI»J=l*7l , 1 = 1 ,71 

63  FORHAT ( IX, 2«au, 7E 16  .71 
03  15  1=1,7 

B9ii>  =a.c 

03  15  <=1.7 

15  B3(ii  = ai (i ) »aoc<,  ii *aa(Ki 

3lN«Ci=SjkT  (B3C1I  • *2  ♦95C2I  • • 2 ♦BB  ( 3 I •*  2 I 

31MR1=SJRT (33  III ••2*93(2I  • • 2 *BB  ( J I ? I 

B2H  = S4RT  (93(41  * * 2*  BBC5)  * • 2 ♦ BB ( 61  • *2  *93(71  **2 1 

NB=1 

IOCT^S 

CALL  L»SQAR(iP,a3,*1,9,M3,'(,K,I3GT,W<AREA,IERI 
WRITE  (5, -,31  (1,33(11  ,1  = 1,31 
4b  F3RJ1AT  ( 1 X,  2H43,  IlO,  Elb.SI 

XS(J0=a9<  11  ♦•203(2I**2»43<3I**2 

X'lRi  = b3RT(  XSQOI 

93L3=X*4R9 

<®=  n 

CALL  jjJbUiJ.Mt 
33  3)  1 = 1,7 

33  A(II  = A(ll*(-ll*M*T  (1 1*^4  7*4(11 

CALL  G 5Jd( 2*<  J» 7, Tl 
03  31  1=1, 7 

31  B(II=3(II-(-ll**I*T  (l  I * . 1 1 3*  B ( I I 
CAL.  GSJd( 3 *<U« 7,  Tl 

33  3?  1=1,7 

32  C(II  = C(I I* ( -1 1 * * I *T  CH*.323*CCI> 

CALL  5i J9< *JQ, 7,  Tl 

03  33  1=1,7 

35  X9(I»=XH<II-(-ll**I*T(II*.J5*X*(ll 
WRITE(b,212T(T(II  , X N( 1 1 ,1  = 1,71 
212  F3R1AT (IX, 1 HT ,E13«S,29XN,£16,SI 
XlA183A  = 51‘,5.  j*3.  93  7 E-!)3 
03  ? <=1,7 
39(0  = XL  ABBDA’XO  Kl 
<<<(  O = 38 ( K I »CC ( < I 
CC(<»  = 93  (< I 
AA(<,  1 1 = A( K I 
AA (A, 2 I =B( Kl 
2 A A ( < , 3 I = C ( < I 
03  SJ  LL  =1  , 3 
X < ( L L I =:  .0 
00  3l  <=1, 7 

63  X<CLLI =XK(LL> *JCCK. LLI *<<< (<l 
X<2=*<(1I**2»X<C2I**2*X<(3I**2 
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RWS  = (XAPSI2*XK2>  7XNRH 
WRITE(b,65)RWS 
65  F3RiAT(lX,4H  RHS,  El  6. 8) 

03  SI  1=1*7 
03  85  J=  1*  3 
GC ( I « J ) = AA  ( I.  Jt 
83  ail . J) =AA< I , J) 

XUTEIb,5D  » 

50  F3RNAT  1 1 HI / IX , ?WA A AN 3 33) 

WRI TE I b , 51 ) < Aft  I < * 1)  « A A ( rC,  2 ) *AA(K«3). 9B  ( <) « ( = 1 * 7) 

51  FSRiATIlX, 3E16.S, E20.S) 

03  11  1=1,7 
0)  11  J= 1*  3 

11  ETII,  J)=QR(I, J) 

03  b 8 1 = 1,7 

03  b 8 3=1*3 
EKI.  J » = C . 0 
GHII, J»=0.0 
03  bS  <=1,3 

£LU.J>=EL(I,J>»ET(I,OMO(K,J) 

68  GNII, J) = GH (I, J) »GGI I, <>  *J3 (K, J) 

03  >3  1=1,7 
03  63  J = 1 , 3 
ETII, Jl=0.0 
GGII, J»=u.O 
03  t>3  <=1,7 

ET(I,J>=£I  II,  JM)UIK,I)*:L(K,  J) 

68  GGU,J)=GG(I.J)»3U(<,  I)*SN(<,3) 

WRITE  lb,  2b)  KET  (I  ,JI,  J=l,  3)  , 1 = 1 ,7) 

2b  F3R1AT(lX,2HET,3Elb.8l 

WRITEtb«27)  IIGSII ,J),J  = 1,3) , 1 = 1 ,7) 

27  F3RNAT I IX, 2HGG,  3E lb ,5) 

33  1?  1=1,7 

03  12  J=l,3 

12  ETII,  J)  = GG(I,  J)-ET(  I,  J) 

WRITE  lb,  2b)  1 1 ET  II  ,J),J  = 1,  3)  ,1  = 1,7) 

33  71  1=1,3 
03  71  J = i,3 

71  Ellll. J)=ET  (I,J) 

03  72  1 = 4,7 

03  72  J=l, 3 

72  E21II-3.JI  =ET  II, J) 

CALL.  LSVALR  I GG, N,  N,  M»  N , 1 3 W , WKAREA,GQ,GU,  3JI 
WRITEIb,29)  ( GQ ( I ) ,1=1,3) 

28  FORMAT I1X, 2M3,«E16, 8) 
3lNRN=ANAXlll.]7GQIl),l.a/GQI2),1.3/G3(3>) 

CALL  LS\/AlR(E1L,N,N,H,N,ISw,W<AREA,G3,GJ,GU 
WRITE  I b, 29 ) IGQII)  ,1=1,3) 

EllNRNr  AMAX1I G3II  ) , GUI  2 ) , GQI 3) ) 

1 = 4 

CALL  LSVALR  IE21,1  ,N,H  ,N  , I SW,W  <AREA  ,GQ,  3U,  Gi/  ) 
WRITE (6,28)  IGQII)  ,1=1,3) 

H = 7 

E21WRN=AMAXl(Ga(l),GQ(2),G3(3)) 

A3C0EF=3INRM*E21NRM 

CU  = A3C3E?/SQRf  II  .3  ♦A9C3£F**2I 


( 


R9S  = SIMRM*t:iiNRM*  I82MM*BINRN*C9  II  2XNRM 

XRMS  = 8INRM*E11MRMMBIMRM**2*B2NM*E21NRM)  / X NRM*  81 NRM  **  2*  E21  NRM*  *2 
*RITE<6,64>R9S.XRHS 

64  FORMAT  ( 1 X,  tH  R9S,  E16.8,U,5M  XRRS.E16.8I 

8RITE  < 6, 22  » RMS,  8INRM,  EilNRM, 3NRM. CHI.XNRM, E21NRM 
22  FORMAT  (IX.  1HM.2E16.  8) 

XKA»SAR=ANRM*BINRM 
KNU=91MRMX  ( ANRI  M OL  Dl 
RRITECS.12) X<A3BAR«XNJ 
12  FORMAT  C1X,  8H<A3!»A8ARtE16. 8, 3H  NU.E16.8I 


33  86  J=l,3 
8(JI=0. 

03  28  1=1,2 

28  H(JI=M(J)»AA(I,  Jl  *AA(I.  Jl 
86  HI  Jl =S3RT(H(JII 

03  83  J= 1. 2 
IllM  = JH 
03  81  I = IL  I M»  3 
SJM  = J . 3 
03  28  <=1,2 

29  SUM=SJN»AA  (K,J)*AA(K,  I> 

80  C33SII,  J)  = SUM/H<  I)  *MIJ)I 

M = 2 
N A A = 3 
N 38  = t 
I AA= I B3=  2 
I 3S  T = 8 

CALL  LL3aAR(AA,B3,H,NAA,RBB,IAA,IBB, I0GT,H<ARE  A, IER) 

MRITElS»25w)(|AA(K»J),J=l»3),K=l»2l 
203  FORMAT  1 1 H j / (IX,  SMPSUEOlNi/,  3E15.2I  I 
HRITE  (6,  lOu  > Hi 

100  FORMATIlHj/lX.SM  IER  =,110) 

MRITE(b.lDl) <<.93<<>,<=1»3) 

101  FORMAT!  IX,  3HJVM,  I 9,  E16. 8) 

SJMt  = 38  < 1) * 8811 1 ♦ 8B 12) *33 (2  > »B8  13 ) *B8<  3) 

XMRM=SOR  T<  SUM1I 
HRl TE IS. 92  > SJMl  , X MR  M 

92  FORMAT  I 180  / 1X,29JVR  S30 , £ 1 6 . 8/1 X , BUS 3R T JWM.Elb.S) 

03  *♦  <=1,2 

44  SI<I=CC(KI-AIIO*3BI1I-8(<)*BB(2I-C<IO*BBI3I 
MRI T £ I 6 < 43  I U.SK  ),  K=l,  2) 

43  FORMAT  I 1HC/  11  X,  5 H RE  SI3.Ild.El6.  8)  ) 

KlHS  = 

L C B 3 ( 1 ) - JB ( 1 1 I **2» IB1I2 I -Q8  < 21  I **2» I BBI 3) -OB<  31  I **2 
MRITEI6.21)  (3B(I»  ,J8II> ,1  = 1,3) 

21  F3R1Ar(lX,*HB818,2E16.8) 

MRI TE I 6, 24  I XLHi.HOLO 
24  FORMAT I1X, 2H**, 2E 16  .81 
XLMS=S1RT(XLHSI f 4 OL  0 
RRITEIS»b2) XlHS.XLRS, RBS.RhS 
62  FORMAT  1 1 X,  *HX  L9  S , El  5.  8 , F 1 5 . 8 , 4*  R MS,  E 1 5. 8,  r 1 5 . 8) 


17 


i 


03  7i  < = it  7 

75  M(K)-S9RT(A«)»AnO  »B<  K » * B ( K ) *C  < X)  *C  ( K) ) 

03  75  1 = 1.6 
JLI1  = IH 
DO  75  J=  JL I M,  7 

C3S=<A(I »*A  <J)»B(  I»  *B(J)*C(I)*C  CJ>)/  ( H < I » * -I  C J»  > 
7b  WRI TE I b , 77 ) I,  J.C3S 

77  FORMAT <1H0 ,5X,tOMROW  COSINE , 213 ,E16.  8) 

00  3L  J=l,2 

1 LI  N=  J *•  1 

03  81  1= IL I M. 3 

81  WRITE (b, 82) J, I, CCOS  (I. J> 

82  FORMAT C 1 HO , 5K , L 3 H COLUMN  COSINE. 2I3.E16.S) 

N=  3 

I A=  I A A 
I SW  = 1 

CALl  LSt/A;.R  U.M.NVI  A.IA.ISW,  MKAREA.S.C.tf) 

WRITE  lb . 88 ) («,$(< ). K=l, 3) 

88  FORMAT ( 1 HO  / (10X.SHSINC  V A L , I 5 , E lb  . 8 ) ) 

WRITE (b« 89)  ( ( G(  I , J)  , J= 1 , N)  .1=1.  M) 

89  F3RMAT  ( 1 HO  / (IX.  HU.  3E15.3)  ) 

WRITE  (b,  93  ) (OMI.  Jl  ,J=1.M)  ,1=1.  N) 

90  FORMAT  ( 1 HO  / (IX.  HV.  3El  5 . 3 ) ) 

WRITE(b,9<»)  ((T<K(I,  J).J  = 1,3)  ,1=1,3) 

94  FORMAT  1 1H0/ (IX, 3-»T<K, 3Elb.  8)  ) 

99(1) =93(2) =39(3)  =0 .0 

CALL  LLSQAR  (TK<,B9,  3,  3 , 1 , 3 , 3 , 8 , WK  ARE  A , I E R > 
WRITE(b*93)  IE R 
93  FORMAT (IX, 5HIER  = ,1  10) 

WRITE<b,95) ((T<<(I, J).J=1,3) ,1=1.3) 

95  FORMAT  (1H0/ (IX,  7-t  TKK  PSI.3E16.8)) 

3L  C3NTINJE 

STOP 

END 


SD3R 3 J T I N£  VG01PI  V,Q,M,N) 

01  HENS  I ON  V(7,7),VA(7,7),U(7,71  , INDEX (13  I , INI/  ( 10)  ,C  (10) 
DIMENSION  WKARE  A ( 25  0) 

IF  IN. GT  . J ) GO  TO  400 
NPO=N»l 
03  53  J = NP Q , M 
03  33  1=1, M 
33  VII, Jl =i.J 
03  2 1=1, H 
03  l J = 1 » N 

1 VAlI.J)=G.o 

2 VIII, I) =1. 3 
03  3 1 = 1 , N 

3 INDEX ( I » - I 
03  51  1=1, N 

IF  IVII,I> • NE • 3 • 1 ) GD  TO  49 
03  S3  J= I, N 

IF  <V< I, Jl .£3.0.3  > GO  TO  30 
I>T0RZ=IN3EX(I) 

INOEXI II =INOEXC Jl 
INOEX I Jl =1  STORE 
GO  T 3 91 
a:  03MUNJE 
WRITEI6, 931 

83  F3R NA  T ( 191/  (IX,  8-INO  PIVOT)) 

STOP 

81  03  92  K=  1,  M 
ST3RE  = V IK, II 
l/K.I)  = V ( < , Jl 
VIX, Jl = 3 TORE 
3T0R£=VA IK, II 
VA  IX,  I>  =Vfl  (<♦  Jl 

82  V 9 I X,  J I = ST  ORE 
49  STORE  = V < I, II 

03  51  J = I , N 
VII,  Jl  =V(I  , JI/3T0RE 
51  VAII,  J)=VA  (I,  JI/’iTORE 
I P3  = I * 1 
03  53  <=  1, H 
IF  K.  Ell.  1 1 GO  ro  5 3 
ST3 R£  = V (K,  1 1 

00  54  J = I » N 

VK,JI=V<<,JI->r3RE*VII,Jl 
54  VAt<,J)=VAt<,JI-3T3RE»VA(I,  Jl 
53  CONTINJE 
5C  OONTINUE 

1 33  T = 8 

CALL  L IN  V?  F I V A , M*  H,  V,  IOGT  ,WKAREA,  IER) 

WRIT£I3,67>  IER 
67  FORMAT (IX, 5 H IER, I 10 1 
03  33  1=1, N 
39  INV IIN3l  XI  I ) ) =1 

HR  I T E I 5 , 2311  ( 1 1 V ( II  ,1  = 1. N) 

201  FORMAT ( 1 X, 3HINV, 1 01 101 

M4ITE  f 6,20OMIN3£  XI  I)  , 1 = 1 , N) 

20u  FORMAT (IX, 5HIN3EX.1 01101 
00  38  J=1,M 


1 


! 


38 
4 0 u 


03  38  1=1, H 
UCI,J»  =y/  (l,INV(J)  ) 
03  4J2  J = M,N 
03  4J1  1=1, M 


401 

402 


an,  j* 


.0 


403 


20 


21 


23 


26 


27 


28 

22 


99 

100 

101 


ac j,i» =i . 

□ 3 43  3 J = 1 , N 
03  433  1=1, N 

a<i.  J>  =v (i , j> 

S JM=  3 . 0 
03  23  1=1, N 
SJM=SUM*Q( 1 , 1 > * ► 2 
SDM=SaRT  (SUM) 

00  21  1=1, M 
3(1,11 =3(1 ,1) />Ji 
03  22  L = 2 , M 
LMO=L-l 
00  23  <= 1, L MO 

cm  =o . o 

03  23  1=1, M 

cm  =o  (<) 4Q  (i > • a(  i, o 

03  26  1=1, M 

03  26  K= 1, L M3 

aci,L)  =j(i,u-;«)*QU,o 

S JM  = 3 . 0 

03  17  1=1, M 

SJM=5J1tQ( I , L ) **Z 

S JM  = S3  R T (SUM) 

03  28  1=1, M 
0(1, L) =0(1 ,L) /SUM 
C3NTIM JE 
03  133  1=1, M 
03  100  J=I,M 
sjM=a. 

03  99  <= 1, M 
SUN=SJM»Q< K,I) *a( K,  J) 

RRI  TE  ( 6,  10  1 ) I , J,  SUM 
F3RM*T(lX,6HIJSJM,2I5,Eli.8> 
RET  JRN 
£90 


ACKNOWLEDGMENTS 


The  author  wishes  to  thank  the  following  individuals  for  their 
interest  and  assistance:  Dr.  Elizabeth  H.  Cuthill  (DTNSRDC,  Code  1805); 
Dr.  Surendra  K.  Dhir  (DTNSRDC,  Code  1844);  and  Mr.  Robert  P.  Eddy 
(DTNSRDC,  Code  1843). 


REFERENCES 

1.  Nobis,  D.  and  C.M.  Vest,  "Statistical  Analyses  of  Error  in 
Holographic  Interferometry,"  Applied  Optics,  vol.  17,  no.  14  (Jul  15,  1978) 
pp.  2198-2204. 

2.  Stewart,  G.W.,  "On  the  Perturbation  of  Pseudo-Inverses, 

Projections  and  Linear  Least  Squares  Problems,"  SIAM  Review,  vol.  19, 
no.  4 (Oct  1977)  pp.  634-862. 

3.  Stewart  G.W.,  "Introduction  to  Matrix  Computations,"  Academic- 
Press,  New  York  (1973). 

4.  Wilkinson,  J.H.,  Rounding  Errors  in  Algebraic  Processes," 
Prentice-Hall,  Englewood  Cliffs,  N.J.  (1961). 

5.  International  Mathematical  & Statistical  Libraries,  Inc., 

"Computer  Subroutine  Libraries  in  Mathematics  and  Statistics,"  1MSL  Library 
3 (CDC  6000  series),  vol.  2,  Houston,  Texas  (Apr  1977). 


INITIAL  DISTRIBUTION 


Copies 

Copies 

Code 

Name 

1 

NRL/Code  8441  (J.  Hansen) 

10 

1844 

D.A.  Gignac 

O 

NAVPGSCOL 

1 

1850 

T.  Corin 

1 59C1/C.  Cantin 

1 M.  Kelleher 

1 

1870 

M.  Zubkoff 

12 

DDC 

2 

1890 

G.R.  Gray 

J.D.  Strickland 

1 

Air  Force  Aero  Res  Lab 

P . Nikolai 

1 

1892 

S . E . Good 

2 

BUS TAN D 

2 

1966 

J.R.  Caspar 

l H.  User 

1 M.  Cordes 

Y.-N.  Liu 

1 

U Maryland/Prof.  G.W.  Stewart 

1 

522.1 

Unclassified  Lib 

(C) 

l 

522.2 

Unclassified  Lib 

(A) 

CENTER  DISTRIBUTION 


Cop ies 

Code 

Name 

1 

1552 

J. 

Bai 

4 

1725 

R.  F 

. Jones,  Jr 

J.E 

. Roderick 

P.N 

. Roth 

L.N 

. G i f f o ra  , . 

1 

1730.5 

J . H 

. Ma 

2 

I 740.4 

M. 

Sal ive 

N.E 

. Fioravant 

1 

1745 

C.  ] 

Ng 

1 

1800 

G . H 

. Gleissner 

1 

1802.2 

F.N 

. Frenkiel 

1 

1805 

E.H 

. Cuthil 1 

2 

1809.3 

D.  1 

Harris 

l 

1820 

A . W 

Camara 

2 

1826 

L.K 

Meals 

L.M 

Culpeper 

1 

1840 

J.W 

Schot 

1 

1843 

R.P 

Eddy 

4 

1844 

S.K 

Dhir 

G.C 

Evers t ine 

F.M 

Henderson 

M.M 

Hurwitz 

DTNSRDC  ISSUES  THREE  TYPES  OF  REPORTS 

1 DTNSRDC  REPORTS.  A FORMAL  SERIES,  CONTAIN  INFORMATION  OF  PERMANENT  TECH 
NICAL  VALUE  THEY  CARRY  A CONSECUTIVE  NUMERICAL  IDENTIFICATION  REGARDLESS  OF 
THEIR  CLASSIFICATION  OR  THE  ORIGINATING  DEPARTMENT 

2 DEPARTMENTAL  RF PORTS.  A SEMIFORMAL  SERIES.  CONTAIN  INFORMATION  OF  A PRELIM 
INARY  TEMPORARY  OR  PROPRIETARY  NATURE  OR  OF  LIMITED  INTEREST  OR  SIGNIFICANCE 
THEY  CARRY  A DEPARTMENTAL  ALPHANUMERICAL  IDENTIFICATION 

3 TECHNICAL  MEMORANDA  AN  INFORMAL  SERIES  CONTAIN  TECHNICAL  DOCUMENTATION 
OF  LIMITED  USE  AND  INTEREST  THEY  ARE  PRIMARILY  WORKING  PAPERS  INTENDED  FOR  IN 
TERNAL  USE  THEY  CARRY  AN  IDENTIFYING  NUMBER  WHICH  INDICATES  THEIR  TYPE  AND  THE 
NUMERICAI  CODE  OT  THE  ORIGINATING  DEPARTMENT  ANY  DISTRIBUTION  OUTSIDE  DTNSRDC 
MUST  BE  APPROVED  BY  THE  HEAD  OF  THE  ORIGINATING  DEPARTMENT  ON  A CASE  BY  CASE 
BASIS 


